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Abstract 

Convex polyhedra are the basis for several abstractions used in static analysis and 
computer-aided verification of complex and sometimes mission critical systems. For 
such applications, the identification of an appropriate complexity-precision trade- 
off is a particularly acute problem, so that the availability of a wide spectrum of 
alternative solutions is mandatory. We survey the range of applications of polyhedral 
computations in this area; give an overview of the different classes of polyhedra 
that may be adopted; outline the main polyhedral operations required by automatic 
analyzers and verifiers; and look at some possible combinations of polyhedra with 
other numerical abstractions that have the potential to improve the precision of 
the analysis. Areas where further theoretical investigations can result in important 
contributions are highlighted. 
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1 Introduction 



The application of polyhedral computations to the analysis and verification 
of computer programs has its origin in a groundbreaking paper by Cousot 
and Halbwachs [43] . There, the authors applied the theory of abstract inter- 
pretation [38,40] to the static determination of linear equality and inequality 
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relations among program variables. In essence, the idea consists in interpret- 
ing a program (as will be explained in more detail in Sections 12.11 and |) on a 
domain of convex polyhedra instead of the concrete domain of (sets of vectors 
of) machine numbers. Each program operation is correctly approximated by a 
corresponding operation on polyhedra and measures are taken to ensure that 
the approximate computation always terminates. At the end of this process, 
the obtained polyhedra encode provably correct linear invariants of the ana- 
lyzed program (i.e., linear equalities and inequalities that are guaranteed to 
hold for each program execution and for each program input). 

As we show in this paper, relational information concerning the data objects 
manipulated by programs or other devices is crucial for a broad range of appli- 
cations in the field of automatic or semi-automatic program manipulation: it 
can be used to prove the absence of certain kinds of errors; it can verify that 
certain processes always terminate or stabilize; it can pinpoint the position 
of errors in the system; and it can enable the application of optimizations. 
Despite this, due to the lack of efficient, robust and publicly available imple- 
mentations of convex polyhedra and of the required operations, the line of 
work begun by Cousot and Halbwachs did not see much development until 
the beginning of the 1990s. Since then, this approach has been increasingly 
adopted and today convex polyhedra are the basis for several abstractions used 
in static analysis and computer-aided verification of complex and sometimes 
mission critical systems. For such applications, the identification of an appro- 
priate complexity-precision trade-off is a particularly acute problem: on the 
one hand, relational information provided by general polyhedra is extremely 
valuable; on the other hand, its high computational cost makes it a fairly 
scarce resource that must be managed with care. This implies, among other 
things, that general polyhedra must be combined with simpler polyhedra in 
order to achieve scalability. As the complexity-precision trade-off varies con- 
siderably between different applications, the availability of a wide spectrum 
of alternative solutions is mandatory. 

In this paper, we survey the range of applications of polyhedral computations 
in the area of the analysis and verification of hardware and software systems: 
we describe in detail one important — and historically, first — application of 
polyhedral computations in the field of formal methods, the linear invariant 
analysis for imperative programs; we provide an account of linear hybrid sys- 
tems that is based directly on polyhedra; and we explain with an example how 
polyhedral approximations can be applied to analog systems. The paper also 
provides an overview of the main polyhedral operations required by these ap- 
plications, brief descriptions of some of the different classes of polyhedra that 
may be adopted, depending on the particular context, and a look at some pos- 
sible combinations of polyhedra with other numerical abstractions that have 
the potential to improve the precision of the analysis. Areas where further the- 
oretical investigations can result in important contributions are highlighted. 
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The plan of the paper is as follows. Section [2] introduces the required no- 
tions and notations, including a minimal exposition of the main concepts of 
abstract interpretation theory. Section [3] demonstrates the use of polyhedral 
computations in the specification of a linear invariant analysis for a simple 
imperative language; a few of the many applications for the analysis of com- 
puter programs are briefly recalled. Section [4] is devoted to polyhedral ap- 
proximation techniques for hybrid systems, which, as shown in Section [5] can 
also be applied to purely analog systems. Section [6] presents several families 
of polyhedral approximations that provide a range of different solutions to 
the complexity/precision trade-off. The most important operations that such 
approximations must provide in order to support analysis and verification 
methods are illustrated in Section [3 Section concludes. 



2 Preliminaries 

We assume some basic knowledge about lattice theory [27|. Let (5, E) and 
(T,d) be two partially ordered sets; the function /: 5 — > T is monotonic 
if, for all x ,xi G S, x E x\ implies f(x ) ^ f(xi). If (5, E) = (T, dt), so 
that /: S — > S, an element x G S such that x = f(x) is a fixpoint of /. 
If (S, E, -L; T, U, n) is a complete lattice, then / is continuous if it preserves 
the least upper bound of all increasing chains, i.e., for all Xq E x\ □ • • • in 
S, it satisfies f[Ux^j = Uf(%i)', m such the least fixpoint of / with 

respect to the partial order denoted lfp/, can be obtained by iterating 
the application of / starting from the bottom element _L, thereby computing 
the upward iteration sequence 



For each fa: So — > T and f\: S± — > T\, the function /o[/i]: (So U Si) 
(T U Ti) is defined, for each x G 5b U Si, by 



For n > 0, we denote by v = (u , . . . ,f n _i) G M™ an n-tuple (vector) of real 
numbers; R + is the set of non-negative real numbers; (v, w) denotes the scalar 
product of vectors v, w G M. n ; the vector OeP has all components equal to 



i- = f°{±) E f\±) E f(i-) E • • • E f(±) E • • • , 



up to the first non-zero limit ordinal tu; namely, 



ifo/ = n-L) = urn)- 
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zero. We write v :: w to denote the tuple concatenation of v 6 1" and w G M. m , 
so that v :: w G R n+m . 

Let x be an n-tuple of distinct variables. Then (3 = f(a, x) txi bj denotes 
a linear inequality constraint, for each vector a G M n , where a / 0, each 
scalar 6 6 R, and ix G {>,>}. A linear inequality constraint (3 defines a 
(topologically closed or open) affine half-space of R n , denoted by con({/3}). 

A set V C R n is a (convex) polyhedron if and only if V can be expressed 
as the intersection of a finite number of affine half-spaces of M n , i.e., as the 
solution V = con(C) of a finite set of linear inequality constraints C (called 
a constraint system). The set of all polyhedra on the vector space M. n is de- 
noted as F n . When partially ordered by set-inclusion, convex polyhedra form 
a lattice (P n , C, 0, M n , l±l, fi) having the empty set and W 1 as the bottom and 
top elements, respectively; the binary meet operation, returning the great- 
est polyhedron smaller than or equal to the two arguments, is easily seen to 
correspond to set-intersection; the binary join operation, returning the least 
polyhedron greater than or equal to the two arguments, is denoted 'l±J' and 
called convex polyhedral hull (poly- hull, for short). In general, the poly- hull of 
two polyhedra is different from their convex hull |110| . 

A relation ip C W 1 x R n (of dimension n) is said to be affine if there exists 
IgN and a,, Cj G M n , h G K and X* G {>, >}, for each % — 1, . . . , £, such that 

i 

Vv, w G R n : (v, w) G tjj ^ A (( c " w ) ( a " v ) + h i) ■ 

i=i 

Any affine relation of dimension n can thus be encoded by I linear inequal- 
ities on a 2n-tuple of distinct variables x::x' (playing the role of v and w, 
respectively), therefore defining a polyhedron in P 2n . The set of polyhedra P n 
is closed under the (direct or inverse) application of affine relations: i.e., for 
each V G P„ and each affine relation f C 1" x 8", the image ip(V) and the 
preimage ■0~ 1 (7- > ) are in P n . 

2.1 Abstract Interpretation 

The semantics of a hardware or software system is a mathematical description 
of all its possible run-time behaviors. Different semantics can be defined for the 
same system, depending on the details being recorded. Abstract interpretation 
[38|39II40| is a formal method for relating these semantics according to their 
level of abstraction, so that questions about the behavior of a system can be 
provided with sound, possibly approximate answers. 

The concrete semantics c G C of a program is usually formalized as the least 
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fixpoint of a continuous semantic function T : C — > C, where the concrete do- 
main (C, C, _L, T, U, n) is a complete lattice of semantic properties; in many 
interesting cases, the computational order 'C' corresponds to the approxima- 
tion relation, so that c% C c 2 holds if Ci is a stronger property than c 2 (i.e., c 2 
correctly approximates c\). 

For instance, the run-time behavior of a program may be defined in terms of 
a transition system (E, t, i), where E is a set of states, t C E is the subset of 
initial states, and t G p(E x E) is a binary transition relation mapping a state 
to its possible successor states. Letting E* denote the set of all finite sequences 
of elements in E, the initial history of a forward computation can be recorded"*! 
as a partial execution trace r = o"o • • • <J m £ E* starting from an initial state 
0o G <- and such that any two consecutive states a, and (Tj + i are related by 
the transition relation, i.e., 0^+1) G £. In such a context, an element of the 
concrete domain fp(E*), C, 0, E*, U, DJ is a set of partial execution traces and 
the concrete semantics is lfp(jF), where the semantic function is defined by 

T = XX G p(E*) .XU{rGE*|r = a Gi} 

U [ ro-j +1 G E* t = cr • • ■ Oi G X, (a*, cr i+ i) G t }. 



An abstract domain_j (D", C, _L, U) can be often modeled as a bounded join- 
semilattice, so that it has a bottom element _L and the least upper bound 
d\ U d\ exists for all d\,d\ G DK This domain is related to the concrete do- 
main by a monotonic and injective concretization function 7: — > C. Mono- 
tonicity and injectivity mean that the abstract partial order is equivalent to 
the approximation relation induced on D* by the concretization function 7. 
Conversely, the concrete domain is related to the abstract one by a partial 
abstraction function a: C >— > D* such that, for each c G C, if a(c) is defined 
then c C 7^a(c)^. In particular, we assume that a(_L) = _L is always defined; 
when needed or useful, we will require a few additional properties. 

For example, a first abstraction of the semantics above, typically adopted 
for the inference of invariance properties of programs [39f40] . approximates 
a set of traces by the set of states occurring in any one of the traces. The 
reachable states are thus characterized by elements of the complete lattice 
(p(E), C, 0, E, U, fl), which plays here the role of the abstract domain. The 
concretization function relating D* = jp(E) to C = p(E*) is defined, for each 
S G p(E), by 

7 (d») = { r G E* I r = a ■ ■ ■ a m , Vz = 0, . . . , m : <7j G d* }. 



1 This is just one of a wide range of possible semantics; by the same approach, other 
semantics may be described and related by abstract interpretation [40] , 

2 To avoid notational burden, we will freely overload the lattice-theoretic symbols 

'U', etc., exploiting context to disambiguate their meaning. 
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The concrete semantic function T: £>(£*) — > p(E*) can thus be approximated 
by the monotonic abstract semantic function A: £>(£) — >■ p(E) defined by 

.A = Ad" G p(E) . $ U 6 U { a' G E | 3a G d tt . (a, a') G t }. 

This abstract semantic function is sound with respect to the concrete semantic 
function in that it satisfies the local correctness requirement 

Vc G C : Vd« 6D« :cC 7(d tt ) ^(c) □ 7(-4(d»)), 

ensuring that each iteration JP(_L) in the concrete fixpoint computation is 
approximated by computing the corresponding abstract iteration A*(a:(_L)). 
In particular, the least fixpoint of T is approximated by any post-fixpoint of A 
[38130] . i.e., any abstract element S G D$ such that A(S) C d tt . 

Actually, the abstraction defined above satisfies an even stronger property, 
in that the abstract semantic function A is the most precise of all the sound 
approximations of T that could be defined on the considered abstract domain. 
This happens because the two domains are related by a Galois connection [39] . 
i.e., there exists a total abstraction function a: C — > D" satisfying 

Vc G C : Vd tt G : ar(c) C d tt c □ 7(d tt ). 



Namely, for all c G p(E*), we can define 



a(c) = f [ (jj G E r = a • ■ ■ cr m G c, i G {0, . . . , m} }. 



For Galois connections it can be shown that a(c) is the best possible approx- 
imation in -D" for the concrete element c G C; similarly, a o JF o 7 (i.e., the 
function ^4 defined above) is the best possible approximation for T [35] . Such 
a result is provided with a quite intuitive reading; in order to approximate 
the concrete function T on an abstract element $ G we first apply the 
concretization function 7 so as to obtain the meaning of d"; then we apply the 
concrete function JF; finally, we abstract the result so as to obtain back an 
element of DK 

Abstract interpretation theory can thus be used to specify (semi-) automatic 
program analysis tools that are correct by design. Of course — due to well- 
known undecidability results — any fully automatic tool can only provide par- 
tial, though safe answers. 
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2.2 Abstract Domains for Numeric and Boolean Values 



The reachable state abstraction described above is just one of the possible 
semantic approximations that can be adopted when specifying an abstract 
semantics. A further, typical approximation concerns the description of the 
states of the transition system. Each state a G £ may be decomposed into, e.g., 
a set of numerical or Boolean variables that are of interest for the application 
at hand; new abstract domains can be defined (and composed [39]) so as to 
soundly describe the possible values of these variables. 

As an expository example that will be also used in the following sections, 
assume that part of a state is characterized by the value of an integer vari- 
able. Then, the domain (p(£), Q, 0, S, U, flj can be abstracted to the concrete 

domain of integers (jp(Int), C, 0, Int, U, PI J . This domain is further approx- 
imated by an abstract domain (int", C, _L, Ll) , via the concretization func- 
tion 7i: Int" — > jp(Int). Elements of Int" are denoted by m", possibly sub- 
scripted. We assume that the partial abstraction function a\\ jp(Int) >— > Int" 
is defined on all singletons {m} G jp(Int) and on the whole set Int. We 
also assume that there are abstract binary operations '©', '0' and '©' on 
Int" that are sound with respect to the corresponding operations on jp(Int) 
which, in turn, are the obvious pointwise extensions of addition, subtrac- 
tion and multiplication over the integers. More formally, for '0' we require 
7i(m5©m») D {m + mi \ m Q G 7 i(™oW G lM) } for each rr&M G Int", 
i.e., soundness with respect to addition. Similar requirements are imposed on 
'©' and '©'. Even though the definition of Int" is completely general, families 
of integer intervals come naturally to mind for this role. 

Suppose now that some other part of the state is characterized by the value 
of a Boolean expression. Then, the domain fp(£), C, 0, E, U, <lj can be ab- 
stracted to the finite domain (p(Bool), C, 0, Bool, U, n) , where Bool = {ff , tt} 
is the set of Boolean values. In general, such a finite domain may be fur- 
ther approximated by an abstract domain (Bool", _L, T, U, n), related to 
the concrete domain by a Galois connection. Elements of Bool" are denoted 
by t", possibly subscripted, and we can define abstract operations '0', '©' 
and '©' on Bool" that are sound with respect to the pointwise extensions 
of Boolean negation, disjunction and conjunction over p(Bool). For instance, 
for the operation '©' to be sound with respect to disjunction on p(Bool), 
it is required that, 7b(4 © A) { V ti to G 7b(4)>*i G 7b(4) } f° r 
each tg and t\ in Bool". Likewise for '©'. For '0' the correctness require- 
ment is that, for each t" in Bool", 7 B (©t") 3 { ->t t G 7 B (t") }• Abstract 

comparison operations ©, © : Int" x Int" — > Bool" can then be defined to cor- 
rectly approximate the equal-to and less-than tests: for each m ,m" G Int", 
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7b("4 © m i) =5 { m o — m i m o G 7i( m o)) m i e 7i( m i) }; likewise for '©'. 

Simple abstract domains such as the ones above can be combined in different 
ways so as to obtain quite accurate approximations |39| . In some cases, how- 
ever, the required precision level may only be obtained by a suitable initial 
choice of the abstract domain. As a notable example, suppose that some part 
of the state a G £ is characterized by n (integer or real valued) numeric vari- 
ables and the application at hand needs some relational information about 
these variables. In such a context, an approximation based on a simple con- 
junctive combination of n copies of the domain Int" described above will be 
almost useless. Rather, a new approximation scheme can be devised by model- 
ing states using the domain (|p(M n ), C, 0, M n , U, , where each vector v G K n 
is meant to describe a possible valuation for the n variables. A further ab- 
straction should map this domain so as to retain some of the relations holding 
between the values of the n variables. If a finite set of linear inequalities pro- 
vides a good enough approximation, then the natural choice is to abstract this 
domain into the abstract domain of convex polyhedra (P„, C, 0, R n , l+l, fl) [43] . 
In this case, the concrete and abstract domains are not related by a Galois 
connection and, hence, a best approximation might not existU] Nonetheless, 
the convex polyhedral hull (partial) abstraction function l±J: p(M. n ) >— » F n is 
defined in most of the cases of interest and provides the best possible approx- 
imation. Most of the arithmetic operations seen before can be encoded (or 
approximated) by computing images of affine relations. 



2.3 Widening Operators 



It should be stressed that, in general, the abstract semantics just described is 
not finitely computable. For instance, both the domain of convex polyhedra 
and the domain of integer intervals have infinite ascending chains, so that the 
limit of a converging fixpoint computation cannot generally be reached in a 
finite number of iterations. 

A finite computation can be enforced by further approximations resulting in 
a Noetherian abstract domain, i.e., a domain where all ascending chains are 
finite. Alternatively, and more generally, it is possible to keep an abstract 
domain with infinite chains, while enforcing that these chains are traversed in 
a finite number of iteration steps [41]. In both cases, termination is usually 
achieved to the detriment of precision, so that an appropriate trade-off should 
be pursued. Widening operators [37|38II40I[4T] provide a simple and general 
characterization for the second option. 



3 This happens, for instance, when approximating an n-dimensional ball with a 
convex polyhedron. 
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Definition 2.1 The partial operator V : £>" x >— » D» «s a widening i/: 

J /or a// c/ tt , e" G^d^e 1 «mp/«es too* $ V e tt «s denned and e i C d tt V e tt ; 
/or a// increasing chains e jZ e| □ ■ ■ • , £/ie increasing chain defined by 
do = ef, and df +1 == d\ V (df U ef +1 ), for i <E ~N, is not strictly increasing. 

It can be proved that, for any monotonic operator A: D* — ► the upward 
iteration sequence with widenings starting at the bottom element d = i_ and 
defined by 



converges to a post-fixpoint of A after a finite number of iterations [4l] . 
Clearly, the choice of the widening has a deep impact on the precision of 
the results obtained. Designing a widening which is appropriate for a given 
application is therefore a difficult (but possibly rewarding) activity. 



3 Analysis and Verification of Computer Programs 

In this section we begin a review of the applications of polyhedral computa- 
tions to analysis and verification problems starting with the the work of Cousot 
and Halbwachs |43I65| . These seminal papers on the automatic inference of 
linear invariants for imperative programs constituted a major leap forward for 
at least two reasons. First, the polyhedral domain proposed by Cousot and 
Halbwachs was considerably more powerful than all the data-flow analyses 
known at that time, including the rather sophisticated one by Karr which was 
limited to linear equalities [80|94] . Secondly, the use of convex polyhedra as an 
abstract domain established abstract interpretation as the right methodology 
for the definition of complex and correct program analyzers. 

We illustrate the basic ideas by partially specifying the analysis of linear in- 
variants for a very simple imperative language. The simplicity of the language 
we have chosen for expository purposes should not mislead the reader: the 
approach is generalizable to any imperative (and, for that matter, functional 
and logic) language [TTJ. The abstract syntax of the language is presented in 
Figure [TJ The basic syntactic categories, corresponding to the sets Int, Bool 
and Var, are defined directly. From these, the categories of arithmetic and 
Boolean expressions and of statements are defined by means of BNF rules. 
Notice the use of syntactic meta- variables: for instance, to save typing we will 
consistently denote by s, possibly subscripted or superscripted, any element 
of Stmt. 

The concrete semantics of programs is formally defined using the natural se- 
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Integers m G Int = Z 
Booleans t G Bool = {tt,ff} 

Variables x G Var = f {xq, x±, X2, ■ ■ ■} 
Arithmetic expressions 



Aexp 3 a ::= m \ x \ a + a± \ ao — a\ \ a * cti 



Boolean expressions 



Bexp 3 b:: 



t I Go — CLi J CLq < CLi 



Statements 

Stmt 3 s ::= skip | x := a | So! s i I if then so else si | while b do s 

Fig. 1. Abstract syntax of the simple imperative language 

mantics approach [79]. This, in turn, is a "big-step" operational semantics 
defined by structural induction on program structures in the style of Plotkin 
[97] . First we define the notion of store, which is any mapping between a finite 
set of variables and elements of Int. Formally, a store is an element of the set 



and denoted by the letter a, possibly subscripted or superscripted. The store 
obtained from a G Store by the assignment of m G Int to x G dom(cr), denoted 
by a[m/x], is defined as follows, for each x' G dom(cr): 



The concrete evaluation relations that complete the definition of the concrete 
semantics for our simple language are defined by structural induction from a 
set of rule schemata. The evaluation relations for terminating computations 
are given by A C (Aexp x Store) x Int, for arithmetic expressions, — > C 
(Bexp x Store) x Bool, for Boolean expressions, and A C (Stmt x Store) x Store, 
for statements. The judgment (a, a) A m means that when expression a is 

executed in store a it results in the integer m. The judgment (6, a) A t 
is similar. Note that expressions do not have, in our simple language, side 
effects. The judgment (s, a) A a' means that the statement s, executed in 
store cr, results in a (possibly modified) store a'. The rule schemata, in the 
form concision , that define these relations are given in Figure El Rule instances 
can be composed in the obvious way to form finite tree structures, representing 
finite computations. Figure [3] shows one such tree. 



Store = { a: V Int | V C Var, V finite } 
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(a ,a)^m (a!, a) 



mi 



(m, cr) A m (x, a) A a(x) (a + ax, cr) A m + mi 

{ao,cr} Am (ai, cr) Ami (a ,ff) Am (ai,cr)Ami 

(a — a\,a) A m — m! (a * a 1; a - ) A m • m 1 

(a ,cr)Am (a 1 ,cr)-^m 1 (a ,cr)Am (a 1 ,cr)-^m 1 

(ao = ax, a) — > (m = mi) (ao < ai, a) — > (m < mi) 



(a, cr) A m 



( So ,Cr}Acr" ( Sl)C r") 



(skip, cr) A cr := a, cr) A a[m/x] 



1,0) -> tt (s ,cr) 



s / 

— > cr 



cr 



(if 6 then s else sx, cr) Act' 



(s ; si,a) -> cr' 

(fe,a)Aff ( Sl ,a)A</ 
(if6 then s else si,cr) Act' 



M Aff 



(while 6 do c, cr) A cr 



(b, Cr) A tt ( C ,Cr)Acr" (while b do C, Cr") A a' 

(while b do c, cr) A cr' 



Fig. 2. Concrete semantics rule schemata for the finite computations of the simple 
imperative language 

The possibly infinite set of all finite trees is obtained by means of a least fix- 
point computation, corresponding to the classical inductive interpretation of 
the rules in Figure El The rule schemata in Figure 0] can be used to directly 
model non-terminating computations and need to be interpreted coinductively 
[42ll86fl04| . The judgment (s, a) A means that the statement s diverges when 
executed in store a. By a suitable adaptation of the computational ordering, 
both sets of finite and infinite trees can be jointly computed in a single least 
fixpoint computation [12.80111 OH . While these semantics characterizations con- 
tain all the information we need to perform a wide range of program reasoning 
tasks, they are generally not computable: we have thus to resort to approxi- 
mation. 

Following the abstract interpretation approach, as instantiated in |102fl03|104j . 
the concrete rule schemata are paired with abstract rule schemata that cor- 
rectly approximate them. Before doing that, we need to formalize abstract 
domains for each concrete domain used by the concrete semantics. 

For simple approximations of integers and Boolean expressions, we consider 
the abstract domains Int" and Bool" introduced in Section I2.2L The last (and 
most interesting) abstraction we need is one that approximates sets of stores. 
We thus require an abstract domain (Store", C, A L-l) that is related, by means 
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(xi, (7 ) A 1 (2, cr ) A 2 (x , <Ti) A 1 (xi, <Ti) A 3 



fa + 2, (T ) A 3 (x -xi,cri) A -2 (0, cr 2 ) A (x ,(T 2 )A-2 



(0,cr )A0 (x ,(T )Al (xi := xi + 2,cr ) A (Ti (x := x -Xi,(Ti) A o- 2 (0 < x , (T 2 ) A ff 



to 

Legend: 



= f {(x , 1), (xi, 1)}, 
o\ = f {(ar ,l),(a:i,3)}, 
^2 d = {(xo,-2),(xi,3)}, 

w = f (while < x do (xi := X\ + 2; x := x — x\ 



(0 < x , (T ) A tt ^(xi := Xi + 2;x := x - £i),<7 ) A a 2 cr 2 ) A a 2 

^while < x do (xi := X\ + 2; x := x — x±), (T ) A <j 2 



Fig. 3. The tree representing a concrete execution of a program 



(sq,<t) 



(s ; si, a) A 
, a) A tt (s , a) 



(s Q ,a)^a' (sx,a') 
i \ 00 

(&,<r)Aff (s u a) 



(if6 then s else si,cr) A (if 6 then s else si, cr) A 

(6, a) A tt (c, cr) A (6, cr) —> tt (c, cr) A cr' (while 6 do c, a') 
(while 6 do c, cr) A (while 6 do c, a) A 



Fig. 4. Additional concrete semantics rule schemata for the infinite computations of 
the simple imperative language 

of a concretization function 73 such that 7s (-L) = 0, to the concrete domain 
(p(Store), C, 0, Store, U, n) . Elements of Store" are denoted by a*, possibly 
subscripted. The abstract store evaluation and update operators 

■[■]: (Store* x Aexp) -> Int*, 
•[• := ■] : (Store* x Var x Aexp) — > Store", 
■[■/■]: (Store* x Var x Int*) -> Store* 



are assumed to be sound with respect to their concrete counterparts, i.e., such 
that, for each cr* G Store*, a G Aexp, x G Var and m* G Int*: 

7i(a*[a]) D j m G Int o G 78(0""), (a, cr) A m }, 
7s (a* x := a]) D j cr' G Store a G 7s(cr*), (x := a, cr) A cr' j, 
7s(cr* m*/x]) D |cr[m/x] G Store cr G 7s(cr*),m G 7i(m") j. 

We also need computable "Boolean filters" to refine the information contained 
in abstract stores, i.e., two functions 4>tt,4>s'- Store* x Bexp — > Store* such 
that, for each cr* G Store* and b G Bexp, 

7s(0 tt (a*,6)) D {a G 7S (<7«) | (6, a) A tt }, 
7s(0ff(^,6)) d{ctG 7 s(^) J (&,cr) Aff}. 

We are now in a position to present, in Figure [5], a possible set of domain- 
independent abstract rule schemata. These schemata allow for the free ap- 
proximation of the 'W right-hand sides in the conclusions. This means that 
if, e.g., 

premise 

(S,(j) A (jf 
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(do, cr") A ml (a^ff") Am} 
(m,a") A ai({m}) (x, a") A cr"[x] (a + ai,cr") Amjffi m" 

/ tt\ a tt / tt\ a tt / tt\ a tt / tt\ a tt 
(ao,ff") -tmj (ai, cr) m\ (a , cr") ~^ rrip (ai,(j") ^mj 

(ao — ai, cr") -A m| © m" (ao * ai, cr") m ?> ® m i 

(ao, ^) & ml (ax, cr") A m { 

(t, a") A a B ({*}) (a = a 1 ,cr») A mj © m} 

(ao,(i'} Amg (ai,ff") Am} 

(a < ai, a") A m © mf (skip, cr") A cr" 

(a, a") A m" (a, cr") A m" 

(sc := a, cr") A <t"[x := a] (sc := a, cr") A a"[m"/a;] 

/ t)\ s tt / tt\ s tt 

(s ,CrS)^Cr{ (flt, (j;) ~> 

(«o; si, 4) A o-| 

(o,g")At" (s ,M^,b))^4 ( Sl ,M°Kb))^4 
(if&thenso else s±, cr") A Oq U erf 

(6,Cr")At" (c,0 tt (cr",6)) A af (while 6 do C, (jf) A (7* 
(while b do C, Cr") A fa (cr», 6) U Crf 

Notes: 

(i) This rule is used if the domain Store" can capture the assignment precisely 
(e.g., when Store" is a domain of convex polyhedra and a is an affine expression). 
Notice that the premise is intentionally not used: its presence is required in order 
to ensure that the abstract tree approximates the concrete tree in its entirety. 

(ii) This rule is used when (i) is not applicable. 

Fig. 5. Abstract semantics rule schemata for the simple imperative language 
is an instance of some rule, then 

premise 

(s,a) ~» o\ 

is also an instance of the same rule for each a\ such that a\ C cr|. Hence the 
schemata in Figure [5] ensure correctness yet leaving complete freedom about 
precision. The ability to give up some precision, as we will see, is crucial in 
order to ensure the (reasonably quick) termination of the analysis. 
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It is possible to prove that, for each (possibly infinite) concrete tree T built 
using the schemata of Figures [2] and [H for each (possibly infinite) abstract tree 
T" built using the schemata of Figure EJ if the concrete tree root is of the form 
(s, a) — > u\ (when the tree is finite) or (s, a) — > (when the tree is infinite) 
and the abstract tree root is of the form (s, cr") -2* a{ with a E 7s (0""), then 
T" correctly approximates T. This means not only that a± G 7s (erf) (when T 
is finite), but also that each node in T is correctly approximated by at least 
one node in TK In other words, the abstract tree correctly approximates the 
entire concrete computation (see [H] for the details). 



It is worth stressing the observation in |104| that, even when disregarding the 
non-terminating concrete computations, the abstract rules still have to be in- 
terpreted coinductively because most of the finite concrete trees can only be 
approximated by infinite abstract trees; for instance, all abstract trees contain- 
ing a while loop are infinite. Since, in general, we cannot effectively compute 
infinite abstract trees, we still do not have a viable analysis technique. The so- 
lution is to restrict ourselves to the class of rational trees, i.e., trees with only 
finitely many subtrees and that, consequently, admit a finite representation. 



The analysis algorithm is sketched in |102j . For expository purposes, we de- 
scribe here a simplified version that, however, is enough to handle the con- 
sidered programming language features. The algorithm works by recursively 
constructing a finite approximation for the (possibly infinite) abstract subtree 
rooted in the current node (initially, the root of the whole tree). The cur- 
rent node n = (^(p,o~D r n ), where r n is a placeholder for the "yet to be 
computed" conclusion, is processed according to the following alternatives: 



(1) If no ancestor of n has p in the label, the node has to be expanded using 
an applicable abstract rule instance. Namely, descendants of the premises 
of the rule are (recursively) processed, one at a time and from left to right. 
When the expansion of all the premises has been completed, including 
the case when the rule has no premise at all, the marker r n is replaced 
by an abstract value computed according to the conclusion of the rule. 

(2) If there exists an ancestor node m = (p, crj^) ~» r m of n labeled by the 
same syntax p and such that <r| C cr^, i.e., if node n is subsumed by 
node m, then the node is not expanded further and the placeholder r n is 
replaced by the least fixpoint of the equation r n = f m (r n ), where f m is 
the expression corresponding to the conclusion of the abstract rule that 
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was used for the expansion of node mQ 
(3) Otherwise, there must be an ancestor node m = (p, crj^) r m of n labeled 
by the same syntax p, but the subsumption condition cr| C cr^ does not 
hold. Then there are two options: 

(a) if the abstract domain Store" is finite, we proceed as in case (pQ); 

(b) if the abstract domain Store" is infinite, to ensure convergence, a 
widening 'V over Store" can be employe and store a\ in node n 
is replaced by u^V (cx^ UaJ). Then, we proceed again as in case (DD). 

The abstract semantics of Figure [5] and the given algorithm for computing 
a rational abstract tree are fully generic in that any choice for the abstract 
domains Int", Bool" and Store" will result into a provably correct analysis 
algorithm. Focusing on numerical domains, the role of Int" can be played 
by any domain of intervals, so that the operations '0', '0' and '©' are the 
standard ones of interval arithmetic [T]; for instance, [m^m^] [m^mj] = 
[mj) + m^mU + m"]. More sophisticated domains, such as modulo intervals 
[95], are able to encode more precise information about the set of integer 
values each variable can take. For Store", a common choice is to abstract from 
the integrality of variables and consider a domain of convex polyhedra which, 
in exchange, allows the tracking of relational information. With reference to 
Figure rule (i) can be applied directly when the arithmetic expression a = 
(a, x) + b is affine; the corresponding polyhedral operation is the computation 
of the image of a polyhedron by a special case of affine relation ^ C K" x 1", 
called single-update affine function: 

(v,w)eip w k = (a, v) + 6A f\ Wi = Vi. 

0<i<n 
i^k 

Another special case, slightly more general than the one above and called 
single-update bounded affine relation, allows among other things to approxi- 
mate nonlinear assignments and to realize rule (ii) . For fixed vectors a, c e R n 
and scalars b,d € R: 

(v,w)eV -<=^ ( a > v ) + b < w k < (c, v) + d A f\ Wi = Vi. 

0<i<n 

Both the rules for the if-then-else and the while constructs require the Boolean 
filters and least upper bound operations: these are realized by means of inter- 

4 As explained in |102jl04| . the computation of such a least fixpoint (in the context 
of a coinductive interpretation of the abstract rules) is justified by the fact that 
here we only need to approximate the conclusions produced by the terminating 
concrete computations, i.e., by the concrete rules of Figure (2J which are interpreted 
inductively. Also note that the divergence rules of Figure [4] have no conclusion at 
all. 

5 If Store" is infinite but Noetherian, we can choose V = f U as a widening. 
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sections (or the addition of individual constraints) and poly-hulls, respectively. 
These, together with the containment test used to detect the reaching of post- 
fixpoints and the widening (see Section [7j) used to ensure termination of the 
analysis algorithm, are all the operations required for the analysis of our sim- 
ple imperative language. More complex languages require other operations: for 
instance, the analysis of languages with command blocks needs to have the 
possibility of embedding polyhedra into a space of higher dimension, reorga- 
nizing the dimensions, and projecting polyhedra on spaces of lower dimension. 
Other operations are needed to accommodate different semantic constructions 
(e.g., affine preimages for backward semantics), to allow for the efficient mod- 
eling of data objects (e.g., summarized dimensions to approximate the values 
of unbounded collections [59]), and to help scalability (e.g., simplifications of 
polyhedra [53]). 

Figure [6] illustrates an abstract computation that, by following the analysis 
algorithm above, approximates the concrete tree in Figure O intervals and 
polyhedra approximate sets of integers and sets of stores, respectively. The 
initial abstract store is given by the polyhedron P = con({x >l,Xi= 1}), 



which approximates all concrete stores a satisfying a(x ) > 1 and cr(xi) = 
1 including the concrete store a in Figure [3l Consider first the lower tree 
in Figure [3l This corresponds to the stage in the computation when all possible 
instances of case (ED) of the algorithm have been applied. In particular, the 
two leftmost subtrees are derived according to the abstract semantics rules in 
Figure [5] by only using case (pQ) of the algorithm. For the rightmost child which 
has still to be expanded, V is a placeholder for its conclusion. It is also noted 
that, in the root of this tree, since Vl = 0, the final result will be the same 
as the value assigned to V . Since the rightmost child, satisfies the conditions 
of case ( l3bl) of the algorithm, the abstract store V\ must undergo a widening 
computation, yielding the abstract store Q . Thus this node has to be replaced 
by (w, Qo) ~^ V. Consider now the upper tree in Figure [3] which has the root 
(w, Qo) ^>Pas above. The two left-most immediate subtrees are derived, as 
in the lower tree, by only using case (00) of the algorithm. The rightmost child 
is initially given Q as a placeholder for its conclusion. Since this node satisfies 
the conditions for case ([2) of the algorithm, it is not expanded further; and the 
value of Q is obtained by finding the least fixpoint solution for the equation 
Q = Ql y Q ; namely, Q f = con({2x + 3x 1 > 5,x < 0}). Thus in the 

conclusion of the root of the upper tree we have V = Qq t±J Q = Qq. Finally, 
the completed abstract tree can be obtained by replacing the rightmost child 
of the lower tree by the upper tree and the placeholder V in the conclusion of 
the root of the lower tree by Qq. 

Based on suitable variations of the simple linear invariant analysis outlined in 
this section (possibly combined with other analyses), many different applica- 
tions have been proposed in the literature. Examples include the absence of 
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(0, Qo) * [0, 0] (x , Qo) ^ T (xi := x x + 2, Q > A Q f (i) (a* : = x - x 1 , Q' ) A Q x (i) 



(0 < x , Qo) ^ T (fa := x, + 2; x := x - x x ), Q ) A Q x (if, Qx) A Q = q£ 



(0,7?o) A[Q,Q] (g ,^ ) A [l,oo] (g; 1 :=a: 1 + 2,^)^^ (i) (x := x - x u V' Q ) A V x (j) 

(0 < x ,P > ^ a B ({tt}) ((x 1 :=xi+2;o;o:=a;o-a; 1 ) J ^)A7' 1 (w, Pi) A p ([3b]) 



(tw, P ) (Po y?) = (la?)=P 

Legend: 



w = f (while < Xo do (xi := Xi + 2; x := Xo — X\) 
V = con({x > l, Xl = 1}), Qo = V V (P Q W Pi) = con({2x + 3 Xl > 5, a* > 1}), 

= M^o, < x ) = P , Qo = 0tt(Qo, < x ) = con({x > 1, x x > 1} 



^ = MVo, < x ) = 0, Q^ = ff (Qo, < x ) = con({2x + 3xx > 5, x < 0} 

V'q = f con({x > l,xi = 3}), Q' = f con({x > l,x± > 3} 



V\ = f con({x > —2, xi = 3}), Qi = f con({x +xi > 1, x x > 3} 



Notes: 

(i) Rule (i) of Figure [5] is used here. 

(j2j) Case (j2j) of the algorithm is applied here. 

3bl) Case (l3bl ) of the algorithm is applied here. 



Fig. 6. Finite approximation of an infinite abstract computation tree 



common run-time arithmetic errors, such as floating-point exceptions, over- 
flows and divisions by zero [28]; the absence of out-of-bounds array index- 
ing [43|112| . as well as other buffer overruns caused by incorrect string manip- 
ulations [49f5T] : the analysis of programs manipulating (possibly unbounded) 
heap-allocated data structures, so as to prove the absence of several kinds of 
pointer errors (e.g., memory leaks) [59.106J; the computation of input/output 
argument size relations in logic programs [25l58|72| ; the detection of potential 
security vulnerabilities in x86 binaries that allow to bypass intrusion detec- 
tion systems [82]; the inference of temporal schedulability constraints that a 
partially specified set of real-time tasks has to satisfy [48]. All of the above 
are examples of safety properties, whereby a computer program is proved to 
be free from some undesired behavior. However, the computation of invariant 
linear relations is also an important, often indispensable step when aiming 
at proving progress properties, such as termination [36i89fl08| . It should be 
also stressed that the same approach, after some minor adaptations, can be 
applied to the analysis of alternative computation paradigms such as, e.g., 
gated data dependence graphs |78| (an intermediate representation for compil- 
ers) and batch workflow networks [111] (a form of Petri net used in workflow 
management). 



4 Analysis and Verification of Hybrid Systems 



Hybrid systems (that is, dynamical systems with both continuous and discrete 
components) are commonly modeled by hybrid automata [4|53f73| . These, of- 
ten highly complex, systems are usually nonlinear (making them computa- 
tionally intractable as they are). However, linear approximations, which al- 
low the use of polyhedral computations for the model checking operations, 
have been used successfully for the verification of useful safety properties 



In this section, we illustrate, by means of examples, how polyhedral compu- 
tations can be used for verifying simple properties of hybrid automata. The 
examples are all instances of linear hybrid systems, a particular class of hybrid 
systems that can be modeled using polyhedra where the continuous behavior 
is specified by linear constraints over the time-derivatives of the variables. 

Definition 4.1 (Linear hybrid automaton.) A linear hybrid automaton 

( of dimension n) is a tuple 

(Loc, Init, Act, Inv, Lab, Trans) 

where the first component Loc is a finite set of locations. The functions 
Init: Loc — > P n; Act: Loc — > P„ and Inv: Loc — > F n define polyhedra. In 
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particular, for each location £ G Loc; Init(£) specifies the set of possible initial 
values the n variables can take if the automaton starts at I; Act(£) specifies the 
possible derivative values of the n variables, so that, if the automaton reaches 
£ with values given by the vector v, then after staying there for a delay of 
t G M.+, the values will be given by a vector v + tw, where w G Act (I); Inv(£) 
specifies the values that an n-vector v may have at £. The fifth and sixth com- 
ponents provide a set of synchronization labels Lab and a labeled set of affine 
transition relations Trans C Loc x Lab x P 2n x Loc, required to hold when 
moving from the source location (the first argument) to the target location 
(the fourth argument). 

Observe that the only differences between this definition of a linear hybrid 
automaton and those in, for example [3|52[70f73j . are presentational; in par- 
ticular, as we have used polyhedra to represent the linear constraints, there 
is no need to provide, as is the case in these other definitions, an explicit 
component of the system consisting of the set of n variables. 

The synchronization labels Lab are required for specifying large systems. Each 
part of the system is specified by a separate automaton, and then parallel com- 
position is employed to combine the components into an automaton for the 
complete system. This ensures that communication between the automata oc- 
curs, via selected input/output variables, between transitions that have the 
same label. Example 14.41 provides a very simple illustration of parallel com- 
position; formal definitions are available in [4|73] and a larger application can 
be found in [93] . 

A linear hybrid automaton can be represented by a directed graph whose 
nodes are the locations and edges are the transitions from the source to the 
target locations. Each node £ is labeled by two sets of constraints defining the 
polyhedra Inv(£) and Act(£). To distinguish these constraints, if, for example 
x is a variable used for the constraints defining Inv(£), x will be used in the 
constraints defining Act(£)0 In the examples, the initial polyhedron lmt(£) is 
assumed to be empty unless there is an arrow to £ (with no source node) labeled 
by the constraint system defining Init(l). Each edge r = (i,a,V,£') G Trans, 
is labeled by a constraint system C defining V and, optionally, by a which is 
only included where it is used for the parallel composition of automata. Since 
V ^ we specify C by using two n-tuples of variables x and x', which 
are interpreted as usual to denote the variables in the source £ and target £' 
locations, respectively. We also adopt some helpful shorthand notation: x++ 
and x — denote x' = x + 1 and x' = x — 1, respectively; also, constraints of 
the form x' = x are omitted. The following examples, taken (with some minor 
modifications) from [4l70] . illustrate the automata. 



The dot notation reflects the fact that these variables denote the derivatives of 
the state variables. 
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Example 4.2 A graphical view of a water-level monitor automaton is given 
in Figure [7l This models a system describing how the water level in a tank 
is controlled by a monitor that senses the water level w and turns a pump on 
and off. When the pump is off, w falls by 2 cm per second; when the pump is 
on, w rises by 1 cm per second. However, there is a delay of 2 seconds from 
the moment the monitor signals the pump to change from on to off or vice 
versa before the switch is actually operated. Initially the automaton is at £q 
with w = 1 and it is required that 1 < w < 12 at all times. Thus the monitor 
must signal the pump to turn on when w = 5 and signal it to turn off when 
w = 10. 

The automaton illustrated in Figure [7] has 2 dimensions with variables w and 
x, where x denotes the time (in seconds) since the previous, most recent, signal 
from the monitor. There are four locations l{ where i = 0, 1, 2, 3. At £q and 
£\ the pump is on, while at £ 2 and £3 the pump is off. At £\ and £3 the monitor 
has signaled a change to the pump switch, but this has not yet been operated. 
Thus we have: 

Init(4) = con({w = 1}), Init(4) = Init(4) = hdt(£ 3 ) = 0, 
lnv(£ ) = con({w < 10}), Inv(4) = Inv(4) = con({x < 2}), 
Inv(4) = con({w > 5}), Act (4) = Act (4) = con({± = w = 1}), 
Act(£ 2 ) = Act(£ 3 ) = con({x = l,w = -2}). 

There are four transitions Tij = (£i,ai,Vi,£j) £ Trans, where % G {0,1,2,3} 
and j + 1 (mod 4); the affine relations are 





= con 


v 1 


= con 


v 2 


= con 


v 3 


= V X . 



({w = 10, x' = 0,w' = w}), 
(i[x — 2, x' — x, w' = w}^j , 
(^{w = 5, x = 0, w' = w}^j , 



Example 4.3 A graphical representation of an automaton for a simplified 
version of the Fischer protocol is given in Figure El This models mutual exclu- 
sion for a system with two processors P\ and P 2 with skewed clocks x% and x 2 , 
respectively. Each processor has a critical section and, at any one moment in 
time, at most one may be in its critical section. This mutual exclusion is en- 
sured by a version of the Fischer protocol which requires that Pi and P 2 share 
a variable k; a process Pi (i = 1, 2) is only able to enter its critical section if 
k = i and Pi may only write to k if k = 0. However, it takes at most a time 
units, as measured by Pi 's clock for Pi to set the value of k to i and it could be 
that the other process Pj may also have started writing j to k. To avoid any 
resulting conflict, the protocol requires that Pi must wait for a further b time 
units, also measured by Pi 's clock, before checking that k = i still holds. The 
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time b is called the delay time. The protocol ensures mutual exclusion only for 
certain values of a and b which depend on the relative rates of x\ and x 2 . Here 
it is assumed that the rate of x 2 is between 0.9 and 1.1 times that of Xi and 
that, for i = 1, 2, the clock x^ is reset to zero at the start of both the write 
process and the delay time for Pi. 

The automaton illustrated in Figure has 5 dimensions with variables a, b, 
X\, x 2 , k. Note that here, a and b are constant for all runs of the automaton 
and this is indicated in the graph by the inclusion of the derivative constraints 
a = 6 = at every location. There are six locations: £q where P\ is idle; at £\ 
where k = and P\ is in the process of writing to k; at £ 2 where k = 1 and Pi 
waits for the delay time of b time units; at £3 where k = 2 since P 2 managed 
to complete writing to k before the delay time of b had expired; at £4 where 
the process P\ is in the critical section; at £5 where P 2 has set k = 2 and the 
mutual exclusion guarantee is violated. All the functions and transitions for 
these locations are as given in Figured 

Example 4.4 A representation of an automaton for a simple task scheduler 
is given in Figure EE This models a scheduler with two classes of tasks A\ and 
A 2 , activated by interrupts I\ and I 2 . Interrupt I\ (resp., I 2 ) occurs at most 
once every 10 (resp., 20) seconds and activates a task in class A\ (resp., A 2 ), 
which takes 4 (resp., 8) seconds to complete. Tasks in A 2 have priority and 
preempt tasks in A\ . It is required that tasks in A 2 never wait. 

The Scheduler automaton given in Figure M is the parallel composition of 
two component automata: Interrupt which models the assumptions about the 
interrupt frequencies; and Task, which models the execution of the tasks. The 
Interrupt automaton, which has a single location 'Intpt', has variables C\ and 
c 2 ; c% (i — 1, 2) measures the time elapsed since interrupt Ij occurred. The 
Task automaton has three locations: 'Idle' when no tasks are running; and 
'Taskl ' and 'Task2' when tasks in classes A\ (resp., A 2 ) are active. It has, for 
each i = 1, 2, variables Xi, which measures the execution time of task i, and 
ki, which counts the number of pending tasks in class task i. 

The combined Scheduler automaton has variables x\, x 2 , k\, k 2 , c\ and c 2 and 
locations which are elements of the Cartesian product of the sets of locations 
for Interrupt and Task. As Interrupt has just one location, each Task loca- 
tion £ is used to denote the corresponding Scheduler location; here, the initial 
lmt(£), derivative Act (£) and invariant Inv(£) polyhedra for the Scheduler are 
the concatenation of the corresponding component polyhedra for the Task and 
Interrupt automata (informally, a concatenation of polyhedra V G P m and 
Q G P„ can be obtained by first embedding V into a vector space of dimen- 
sion n + m and then add a suitably renamed- apart version of the constraints 
defining Q). Each transition (£ } a,V } £') in the Task automaton not triggered 
by interrupts I\ and I 2 has a transition (£, a, Q, £') in the product automaton 
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w = 10, x' = 
signal pump off 



switch on 



switch off 



to = 5, x' = 



signal pump on 




Fig. 7. Water-level monitor 
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Fig. 8. Fischer protocol (simplified) 



ci > o, co > o 



Interrupt 



h;ci> 10, 



true 
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; c 2 > 20, c' = 



ii;fcHH 



xi = X2 = k\ = ki = 



Task 



x\ = 4, k\ < 1, 
ki — ,x' x = 0/ 



X2 = 8, k2 < 1, fci = 0, 
k2--,x' 2 = 




x\ = 4,fci > 2, fci--, ^ = 



X2 = 8, k2 < 1, fci > 1, 

k2--,x' = 

Taskl " Task2V 



X2 = 8, k2 > 1, k2 , x' 2 = 



Fig. 9. Scheduler 
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where Q £ P6 is obtained by embedding V into a vector space of dimension 
6. Letting i = 1,2, for transitions (£,Ii,V,£') and (Intpt, Ij, V', Intpt) in the 
Task and Interrupt automata, respectively, there is a transition (£,I{, Q,£') in 
the product automaton where Q EF e is obtained by concatenating V and V . 

Given a linear hybrid automaton, the aim of an analyzer is to check, or even 
find sufficient conditions that ensure, that a valid run of the system cannot 
reach a location and vector of values that violate some requirement of the 
system. For instance, in Example 14.21 we need to show that the water level 
always lies between 1 cm and 12 cm; in Example l4.3l we need to find conditions 
on a and b so that at most one processor can be in its critical section at any 
one time; in Example 14.41 we need to show that no task in A 2 will ever wait. To 
show how polyhedral computations can be used to prove such properties, we 
first define more formally such a run and how reachable sets may be computed. 
Note that these definitions follow, with only minor changes, the approach in 

EDI. 



Letting 7i = (Loc, Init, Act, Inv, Lab, Trans) be a linear hybrid automaton in 
n dimensions, a state s of 7i consists of a pair (£, v), where £ £ Loc and 
v £ Inv(f). Given states s = (£,v) and s' = (£', v'), a time delay t £ R+ and 
a vector w £ Act(£), 



>* s 



is a step of Ti provided that, for all t' £ [0, t), v + t'w £ Inv(£) and, for some 
(£,a,V,£') £ Trans, (v + tw) :: v' £ P. A run of TC is a sequence (finite or 
infinite) of steps 

S0^* Sl-^ «2- - (4.1) 

where the initial state so = (^o, v o) satisfies the condition v £ Initio)- An 
infinite run diverges if the sum J2i>o U diverges. For each divergent run given 
by ( 14. ip where, for i > 0, Si = (£i, Vj), we associate a (state) behavior (3 which 
is a total function from time to states: that is, /3(0) = s and, for each t > 0, 
(3(t) = f (£i,v), where 



min< k £ N 



t j > t \ and v = Vj + w» 1 1 - ^ 

i=0 J V j<% 



A state s is reachable if there exists a divergent run with behavior (3 and time 
t £ M + such that (3{t) = s. The set of all reachable values Ri for a location £ 
is defined as: 

R e = { v £ W l | 3t £ R + . pit) = (£, v) }. 

The set of reachable values Ri at a location £ can be characterized by a system 
of fixpoint equations that are defined in terms of sets of reachable values R^ 
at locations £' where (£',a,V,£j £ Trans. These equations use the following 
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operations on sets of vectors in M n . Let V, Q G P2 n and SCI". Then 
fo(S) = { v' G K rt vgS,v::v'gP}; 

S/g d ={v + iwGl n |ve5,weS,teM + }. 

Note that, if S G P„, then also ip v {S) G P n and S / Q G P n . The 
operator, called the t«me elapse operator, was first proposed in [70]. We can 
now provide the fixpoint equation for Rf. 

Rt = ( (lmt{£) U [J ^(-^0 n Inv{£) \ / Act(£) J n Inv(£). (4.2) 

Informally, the fixpoint equation for i?£ says that the reachable values at the 
location £ are obtained by letting the time elapse either from an initial value 
for £ or from a value obtained from an incoming transition. However, the 
fixpoint Equation (14.21) cannot handle strict constraints correctly and needs 
modifying; this is illustrated in the following example. 

Example 4.5 Consider again Example l4.2l Then, just applying Equation (14.2ft 

(as proposed in [69.70],), the sets of reachable values at locations £±,£2,^3 are 
empty. The reason for this is that, for example, at location £0, the strict con- 
straint w < 10 must hold, while in the transition from £q to £\, the transition 
condition w — 10 has to hold. On the other hand, it follows from the definition 
of a step, that since one of the derivative constraints at £q is w = 1; the water 
level w may continue to increase up to the topological closure of Ri which is 
consistent with w — 10. 

To resolve this problem, in Equation ( 14.21) defining the concrete computation, 
Re needs to be replaced by 

c{Re)n(R e / kct{£')), (4.3) 

where c(R' e ) denotes the topological closure of R' e C M. n . 

Observe that, although the linear hybrid automata are specified by means of 
polyhedra, the reachable set Ri for a linear hybrid automaton and location £ 
may not be in the form of a convex polyhedron. Thus, to verify that some states 
of an automaton are unreachable using the standard polyhedral computations, 
approximations are needed. In particular, in the fixpoint Equation ( 14.21) (or 
( 14.31) ), the set operations have to be replaced by the corresponding polyhedral 
operations. In fact all the operations in (14.21) except set union can be used 
as they are, since they transform polyhedra to polyhedra. Just the set union 
operation has to be replaced by the poly- hull operation 'W described in Sec- 
tion El Thus the following fixpoint equation computes an approximation R\ 
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to the reachability set Re. 
R\={ (lnit(£) W |+J ^v(Rl') H lnv(£) \ / Act(£) ) n Inv(£). (4.4) 

As for the concrete fixpoint equation, to correctly handle the strict constraints 
Equation (14.41 ) needs to be modified by replacing R\, with 

c(i2jj,) n (r\, / Act(O)- 



If we let R" denote the tuple { -Rjj | £ G Loc } we can write Equation (14. as 

R\ = F,(B}). 

For all £ G Loc, we write Rf (0) = and, for all k > 1, Rf fc+1) = i^Rf ). 
Then R 11 can be computed iteratively provided the sequence R^°) , R^ 1 ) , . . . 
does not diverge. To handle diverging sequences, we apply a widening (see 
Section I7T21 . Note that we do not have to apply it at all locations. Let W be a 
set of locations that cut all cyclic paths in the graph of the hybrid automaton 
(that is, each loop of the directed graph contains at least one location in W). 
Then the following set of fixpoint equations is guaranteed to converge: 

<4VF t{ R>), iUeW 
e \F e (Rt), if£GLoc\W. V ; 

Example 4.6 Consider again Example 14.21 As there is a single loop passing 
through Eq, it is sufficient to define the set of widening locations asW = {^o}- 

With the modified form of Equation 114.41 ) and the polyhedra widening of [43j ; 
the computation requires three iterations resulting in polyhedra defined by con- 
straint systems C{ for < % < 3 where: 

C = {1 < w < 10}, Ci = {w - x = 10, 10 < w < 12}, 

C 2 = {w + 2x = 16, 5 < w < 12}, C 3 = {w + 2x = 5, 1 < w < 5}. 

Example 4.7 Consider again Example 14.31 The analysis terminates without 
widening in just two iterations with the resulting polyhedron at t§ defined by 
the constraint system: 

C = {k = 2, 10a > 96, < b < x u 9x x < 10^ < Hxi, 

llxi + 10a > 10x 2 + 116}- 

It therefore follows that, to ensure that there can be no run with a state at 
location £5, it is sufficient that 10a < 96. 
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Example 4.8 Consider again Example 14.41 By applying the above mentioned 
polyhedra widening at location 'Task2' only, the analysis for the product au- 
tomaton terminates in four iterations. After projecting away the variables C\ 
and c 2 , the reachable values are given by polyhedra defined by constraint sys- 
tems Cto, Cti, and Ct 2 for locations 'Idle', 'Taskl' and 'Task2', respectively, 
where: 

Cto = {xi = x 2 = h = k 2 = 0}, 

Cti = {0<x l <A J x 2 = 0, k x = 1, k 2 = 0}, 

C t 2 = {x 2 >0,x 2 < 8, Aki > xi, xi >0,k 2 = 1}. 

So it can be concluded that, at each location of the automaton, k 2 < 1 and, 
hence, no task in class A 2 will ever have to wait. However, as noted in [70] . 
because of the convex hull approximation, with the polyhedral domain the ana- 
lyzer fails to show that k\ < 2. We therefore redid the analysis using a domain 
of powersets of polyhedra (see Section \6.2\) and, after taking the poly-hull of 
the final sets and projecting away the variables C\ and c 2 , we obtained the poly- 
hedra defined by constraint systems C' t0 , C' tl andC' t2 for locations 'Idle', 'Taskl' 
and 'Task2', respectively, where: 

C' t0 = {xt = x 2 = k x = k 2 = 0}, 

C' tl = {0 < x x < 4, x 2 = 0, h = l,k 2 = 0}, 

C' t2 = {0 < xi < 4, < x 2 < 8, 4ki >xi> 2ki - 2, 

x x + x 2 > lOJfei - 10, kx < 2, k 2 = 1}. 

This verifies that kx < 2 and k < 1 in every state of any run of the automata. 

Hybrid systems with affine or nonlinear dynamics do not fit the above spec- 
ification of a linear system so that the verification techniques described here 
are not directly applicable. Nonetheless, by partitioning the continuous state 
space and over-approximating the dynamics in each of the partitions, the same 
techniques used to verify linear hybrid automata can be used in these more 
general cases [50|53f74[76[101| . Such an approach has also been successfully 
applied in the verification of analog circuits, as discussed in the following sec- 
tion. 



5 Analysis and Verification of Analog Systems 

The idea of applying formal methods, that originated in the digital world, to 
analog systems was put forward in [71]. This is an important step forward 
with respect to more traditional methods for the validation of analog circuit 
designs. A formal verification tool can, for example, ensure that a design sat- 
isfies certain properties for entire sets of initial states and continuous ranges 
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(a) Circuit schematic 



V d 

(b) Tunnel diode characteristic 



Fig. 10. Tunnel-diode oscillator circuit 
of circuit parameters, something that cannot be done with simulation. 

In [44] and [64], polyhedral approximations were successfully used in the veri- 
fication of analog circuits. Here, we use a simple example, taken from [55], on 
the verification of an oscillator circuit to illustrate the approach To verify 
properties of the (cyclic) behavior of such circuits, cyclic invariants have to be 
determined. To establish a cyclic invariant for a given set of initial states and 
ranges for the circuit parameters, one has to show that the circuit returns to 
a subset of those initial states, which implies the system will keep traversing 
the same states indefinitely. From such an invariant, a number of properties 
of the oscillator can be established |56|. 



Consider the tunnel-diode oscillator schematized in Figure 10(a) The state of 
the system at a given instant of time is completely characterized by the values 
of the inductor current II and the diode voltage drop V d . With these as the 
state variables, the system is described by the second-order state equations 



V d = l/C(-I d (V d )+I L ), (5.1) 
i L = l/L(-V d -RI L + V in ). (5.2) 

In [55] it is shown how a cyclic invariant can be obtained for this circuit using 
the PHAVer system. First, a piecewise affine envelope is constructed for the 



tunnel diode characteristic I d (V d ) depicted in Figure 10(b): for the particular 
example analyzed in |55| . sufficient precision is obtained by dividing the range 
V d G [— 0.1 V, 0.6 V] into 64 intervals, resulting in a piecewise affine model 
of (15.11) . Forward reachability computation with PHAVer can obtain the set 
of states depicted in Figure [TU These are the states reachable from the set of 
initial states corresponding to V d G [0.42 V, 0.52 V] and II = 0.6 mA (the base 
of the downward- facing triangular shape in Figure [HI). Taking into account 
that the loop shape constituted by the reachable states is traversed clockwise, 
it can be seen that the inductor current returns to the initial value of 0.6 mA 
with a diode voltage drop that is well within the initial range [0.42 V, 0.52 V]. 
The set of reachable states so obtained is thus an invariant of the circuit. 



7 For a more general view, we refer the interested reader to the cited literature and 

to EE]. 
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Fig. 11. Reachable states of the tunnel-diode oscillator (dashed) 

In [55] it is shown that, due to over-approximation, forward reachability can 
fail to determine invariants of more complex circuits. A new technique com- 
bining forward and backward reachability with iterative refinement of the 
partitions is thus proposed and shown to be more powerful and efficient. 



6 Families of Polyhedral Approximations for Analysis and Verifi- 
cation 

For several applications of static analysis and verification, an approximation 
based on the domain of convex polyhedra can be regarded as the most appro- 
priate choice. In this section we discuss alternative options (simplifications, 
generalizations, and combinations with other numerical domains) that might 
be considered when trying either to reduce the cost of the analysis, or to 
increase the precision of the computed results. 

6. 1 Simplifications of Polyhedra 

There are contexts where approximations based on the domain of convex poly- 
hedra, no matter which implementation is adopted, incur an unacceptable 
computational cost. In such cases, the static analysis may resort to further 
simplifications so as to obtain useful results within reasonable time and space 
bounds. 

A first, almost traditional approach is based on the identification of suitable 



29 



syntactic subclasses of polyhedra. The abstract domain of bounding boxes (or 
intervals [37]) is based on polyhedra that can be represented as finite con- 
junctions of constraints of the form ±£j < d or ±£j < d, leading to the 
specification of operations whose worst-case complexity is linear in the num- 
ber of space dimensions. As a more precise alternative, the class of potential 
constraints [2 8.23 45 47 83j , also known as bounded differences, allows for con- 
straints of the form Xi — Xj < d or ±Xi < d; the generalization proposed in [21~] . 
also admits constraints of the form Xi + Xj < d, leading to the abstract domain 
of octagons [90J. In these last two cases, the operators are characterized by 
a worst-case time complexity which is cubic in the number of space dimen- 
sions. For all of the approximations mentioned above, improved efficiency also 
follows from the fact that the corresponding computations are simple enough 
to allow for the adoption of floating-point data types: in contrast, the speci- 
fication of safe and efficient floating-point operations for general polyhedra is 
an open problem, so that polyhedra libraries have to be based on unbounded 
precision data types. 

Several alternative (syntactic and/or semantic) simplification schemes have 
been put forward in the recent literature. The Two Variables per Linear In- 
equality abstract domain is proposed in |107j . where constraints take the syn- 
tactic form axi + bxj < d. In jlOQj . an arbitrary family of polyhedra is chosen 
before starting the analysis by fixing the slopes of a finite number of linear 
inequalities, which are called the template constraints; linear programming 
techniques are then used to compute precise approximations in the consid- 
ered class of shapes. In contrast, in [99], general polyhedra are allowed, but 
the corresponding operations (in particular, the poly-hull and the image of 
affine relations) are approximated by less precise variants so as to ensure a 
polynomial worst-case complexity in the size of the inputs. An even more 
flexible approach is proposed in [53], where arbitrary polyhedra are approxi- 
mated, when they become too complex, by limiting the number of constraints 
in their description and/or the magnitude of the coefficients occurring in the 
constraints. These more dynamic approximation schemes are promising, in 
particular for those applications where nothing is known in advance about the 
syntactic form of the constraints that will be computed during the analysis. 

An important observation to be made is that there is no actual need to prefer 
a priori (and therefore commit to) a specific abstract domain: the analysis 
tool may be based on several abstractions, safely switching from more precise, 
possibly costly domains to more efficient, possibly imprecise ones, and vice 
versa, depending on the context. When replacing a generic polyhedron by a 
simpler one, the problem of the identification of a good over-approximation has 
to be solved. Depending on the context, the approaches may vary significantly. 
At one extreme, when efficiency is really critical, the adoption of syntactic 
techniques should be pursued: for an interesting example, we refer the reader 
to one of the simplification heuristics used in [53], where the efficient selection 
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of a small number of linear inequalities out of a constraint system is driven 
by a simple, yet effective reasoning on the measure of the angles formed by 
the corresponding half-spaces. At the other extreme, linear programming (LP) 
optimization techniques may be used so as to obtain the best match in the 
considered class of geometric shapes. For instance, the precise approximation 
of a polyhedron by a bounding box (resp., a bounded difference or octagon) 
can be implemented by a linear (resp., quadratic) number of optimizations of 
a class of LP problems, where the objective function varies while the feasible 
region is invariant and defined by the constraints of the polyhedron. Note 
that, if correctness has to be preserved, it is essential that no rounding error 
is made on the wrong side, so that classical floating-point implementations 
of LP solvers have to be considered unsafe, unless the computed results can 
be certified by some other tool. Alternatively, it is possible to consider LP 
implementations based on unbounded precision data types. 

When the number of space dimensions to be modeled is beyond a given thresh- 
old, the whole analysis space can be split into a finite number of smaller, more 
manageable components, thereby realizing a further simplification scheme that 
can be combined with those described above. The splitting strategy varies con- 
siderably. In [67|68] . Cartesian factoring techniques are used so as to dynami- 
cally partition the space dimensions of a polyhedron into independent subsets; 
the orthogonal factors are then approximated by lower dimensional polyhedra 
with no precision penalty. In an alternative approach described in [28], many 
(possibly overlapping) small subsets of space dimensions, called variable packs, 
are identified before the start of the analysis by means of syntactic conditions; 
the relations holding between the variables in each pack are then approximated 
by using an octagonal abstraction. A variation of this is described in |112| , 
where non-overlapping variable packs are dynamically computed (and possibly 
merged) during the analysis, whereas the relations between the variables in a 
pack are approximated by means of potential constraints. In |112| it is also 
observed that, since the average size of variables packs is small (5 variables), 
more precise approximations based on general polyhedra should be feasible. 



6.2 Generalizations of Polyhedra 

There are applications where the restriction to the domain of convex polyhedra 
is intrinsically inadequate. This may happen, not only when the verification 
property of interest is itself non-convex, but also when the adopted computa- 
tion strategy requires that a convex property is proved by passing through a 
non-convex intermediate approximation. This was the case in Example 14.81 of 
Section [H where the upper bound {k\ < 2) on the number of waiting processes 
for class A x was obtained by switching from the domain of convex polyhedra 
to the domain of finite sets of polyhedra. 
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The finite powerset domain construction [9] is a special case of disjunctive 
completion |39| . a systematic technique to derive an enhanced abstract domain 
starting from an existing one. A finite powerset domain implements disjunc- 
tions by maintaining an explicit (hence finite) and non-redundant collection 
of elements of the base-level domain: non-redundancy means that a collection 
is made of maximal elements with respect to the approximation ordering, so 
that no element subsumes another element in the collection. 

For a better understanding of the concepts, which are described in completely 
general terms in [T7], let us consider the application of the finite powerset 
construction to the domain of convex polyhedra. This instantiation (which 
is the one also adopted for the examples developed in [17]) can be used to 
model nonlinear systems as described, e.g., in Section [5l Then, an element 
of the abstract domain is a finite set of maximal convex polyhedra, so that 
no polyhedron in the set is contained in another polyhedron in the set. The 
powerset domain is a lattice: the bottom and top elements are and {M n }, 
respectively; the meet is obtained by removing redundancies from the set of all 
possible binary intersections of an element in the first powerset with an element 
in the second powerset; while the binary join is the non-redundant subset of 
the union of the two arguments. Most of the other abstract operations needed 
for a static analysis using the finite powerset domain are easily obtained by 
"lifting" the corresponding operations defined on the base-level domain, and 
then reinforcing non-redundancy. For instance, the computation of the image 
of a finite powerset under an affine relation is obtained by computing the 
image of each polyhedron in the collection. However, the construction of a 
provably correct widening operator has only recently been addressed in [17] 
(see Section E2J) . The generic specification of the abstract operators of the 
finite powerset domain in terms of abstract operations on the (arbitrary) base- 
level domain allows for the development of a single implementation which is 
shared by all the possible instances of the domain construction. 

An alternative abstraction scheme has been proposed in [20] for the case of 
finite conjunctions of polynomial inequalities. Intuitively, a polynomial con- 
straint can be approximated by means of a linear constraint in a higher di- 
mension vector space, so that the different terms of the polynomial (e.g., xo, 
x xi, Xq) are mapped to different and independent space dimensions; these 
linear constraints are then used to perform an almost classical linear relation 
analysis based on convex polyhedra. Due to the linearization step, most of the 
precision of the polynomial constraints is initially lost; however, some of the 
relations holding between the different terms of the original polynomial can be 
recovered by adding further constraints that are redundant when interpreted 
in the polynomial world, but do contribute to precision in the linearized space. 
In particular, in [20] the polynomial constraints are mapped into finitely gen- 
erated polynomial cones and a degree-bounded product closure operator is 
systematically applied so as to improve accuracy. As a trivial example, let 
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the polynomial terms Xq, x\ and x§x\ be mapped to the space dimensions yo, 
yi and y2, respectively. Then, the linearization of the polynomial constraints 
x > and X\ > will produce a polyhedron that, while satisfying y > and 
yi > 0, leaves variable 7/2 totally unconstrained. By applying the product clo- 
sure operator we also obtain the linear constraint yi > 0, thereby recovering 
the non-negativity of term XqX\. 



6.3 Combinations with other Numerical Abstractions 

We observe that there are two basic kinds of numerical abstractions for approx- 
imating the values of the program variables: outer limits (or bounds within 
which the values must lie) and the pattern of distribution of these values. 
The first can be approximated by (constructions based on) convex polyhedra, 
while the second can be approximated by sets of congruences defining lattices 
of points we call grids [T0|60|62| . Before considering how these and similar 
domains may be combined, we give a brief overview of the domain of grids. 

Any vector that satisfies (a, v) = b + fif, for some /i G Z, is said to satisfy 
the congruence relation (a, v) =f b. A congruence system K is a finite set of 
congruence relations in M n . A grid is the set of all vectors in M. n that satisfy the 
congruences in fC. The domain of grids G n is the set of all grids in M. n ordered 
by the set inclusion relation, so that the empty set and M n are the bottom 
and top elements of G n respectively and the intersection of two grids is itself a 
grid. Thus, as for the domain of polyhedra, the domain of grids forms a lattice 
(G n , C, 0, R n , l±l, n) where tfcl denotes the join operation returning the least grid 
greater than or equal to the two arguments. For more details concerning all 
aspects of the domain of grids, see [TO] . 

The distribution information captured by grids has a number of applications 
in its own right: for instance, to ensure that external memory accesses obey 
the alignment restriction imposed by the host architecture, and to enable 
several transformations for efficient parallel execution as well as optimizations 
that enhance cache behavior. However, here we are primarily concerned with 
applications that can benefit from the combination of the domain of grids with 
that of convex polyhedra. For instance, knowing the frequency (and position) 
of the points in a grid, we can shrink the polyhedra so that the bounding 
hyperplanes pass through the grid values; if this leads to a polyhedron with 
reduced dimension (such as a single point) or one that is empty, it can lead, 
not only to improved precision, but also a more efficient use of resources by 
the analyzer [5ll96f98] . 

Generic constructions, such as direct and reduced product, can be used to 
provide a formal basis for the combination of the grid and polyhedral domains 
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[39] although the exact choice of product construction used to build the grid- 
polyhedral domain needs further study. Both the direct and reduced products 
have problems: the direct product has no provision for communication between 
the component domains, thereby losing precision; while the reduced product, 
which is the most precise refinement of the direct product, has exponential 
complexity. It is expected that, for grid-polyhedra, the most useful product 
construction will lie between these extremes. For instance, as equalities are 
common entities for both constraint and congruence systems, if an equality 
is found to hold in one component, it is safe to just add this to the other 
component. In addition, in an element of the grid-polyhedral domain, any 
hyperplane that bounds the polyhedron component could be moved inwards 
until it intersects with points of the grid with only linear cost on the number of 
dimensions. Of course, this reduction on its own is not optimal since the grid 
points in the intersection may not lie in the polyhedron itself. For optimality 
or, more generally, so as to gain additional precision, we need to experiment 
with various forms of the branch-and-bound and cutting-plane algorithms [8l] 
already well-researched for integer linear programming. What is needed is a 
range of options for the product construction allowing the user to decide on 
the complexity/precision trade-off. Further work on this is needed, including 
an investigation of other proposals for generic products that lie between the 
direct and reduced product, such as the local decreasing iteration method [6T] 
and the open product construction [34] . 



7 Polyhedral Computations Peculiar to Analysis and Verification 

As observed in the previous sections, the analysis of the run-time behavior 
of a system can be broken down into a set of basic operations on the cho- 
sen abstract domains. This means that each abstract domain should provide 
adequate computational support for such a set and, where appropriate, fur- 
ther operations that might be useful for tuning the cost/precision ratio. In 
this section, we discuss several key issues relevant to the design and imple- 
mentation of an abstract domain of, or based on, convex polyhedra. Before 
going into further detail, it should be stressed that the particular context of 
the application plays a significant and non-trivial role here. For instance, in 
many computational complexity studies, it is assumed that a small number of 
operations (often, just a single one) can have arbitrarily large operands; also, 
it is typically required that exact results have to be computed. These assump- 
tions taken together may be inappropriate in the context of static analysis: it 
is quite often the case that a large number of operations will have only small 
or medium sized operands; also, whenever facing an efficiency issue, the ex- 
actness requirement can be dropped (provided soundness is maintained). As 
a consequence, the evaluation of alternative algorithmic strategies should be 



34 



largely based on practical experimentation. 



7.1 The Double Description Method 

Convex polyhedra are typically specified by a finite system of linear inequality 
constraints and for this representation there are known algorithms (e.g., based 
on Fourier-Motzkin elimination [84|105j ) for most of the operations already 
mentioned. 

An alternative approach is based on the double description method due to 
Motzkin et al. [92]. This method was originally defined on the set of topologi- 
cally closed convex polyhedra, a sub-lattice (CP n , C, 0, R n , l+l, U) of the lattice 
of (not necessarily closed, or NNC) polyhedra P n . In the double description 
method, a closed polyhedron may be described by using a system of non-strict 
linear inequalities or by using a generator system that records its key geomet- 
ric features. The following is the main theoretical result, which is a simple 
consequence of well-known theorems by Minkowski and Weyl |110| . 

Theorem 7.1 The set V C R n is a topologically closed convex polyhedron if 
and only if there exist finite sets i?,PC R™ of cardinality r and p, respectively, 
such that ^ R and V can be generated from (R, P) as follows: 

V = { Rp + Ptt G R™ | p G R;, tt G RP , n=i 7Ti = 1 }. 

Intuitively, a point of a polyhedron V is obtained by adding a convex combi- 
nation of the vectors in P (the generating points) to a conic combination of 
the vectors in R (the generating rays). 

It turns out that constraint and generator descriptions are duals: each repre- 
sentation can be computed starting from the other one. Clever implemen- 
tations of this conversion procedure, improving on the Chernikova's algo- 
rithms [30|31f32j . are the starting point for the development of software li- 
braries that, while being characterized by a worst case computational cost 
which is exponential in the size of the input, turn out to be practically useful. 
A common characteristic of these implementations is the exploitation of in- 
crementality, whereby most of the computational work done for an operation 
is reused to efficiently compute small variations of the corresponding result. 
Further computational enhancements are obtained by the adoption of suitable 
heuristics, ranging from the efficient handling of adjacency information [85j, 
to a careful choice of ordering strategies for the computation of intermediate 
results pi7ll57] ; the overall construction typically relies on a tight integration 
of the basic algorithms with a carefully chosen set of data structures [T6] . 
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An important motivation for the adoption of an implementation based on 
the double description method is that the ability to switch from a constraint 
description to a generator description, or vice versa, can be usefully exploited 
to provide simple implementations for the basic operations on polyhedra. For 
instance, set intersection is easily implemented by taking the union of the 
constraint systems representing the two arguments, whereas the poly-hull is 
implemented by joining the generator systems representing the two arguments; 
and the test for emptiness can be implemented by checking that the generator 
system has no points. Moreover, a test for subset inclusion V C Q can be 
implemented by checking if each point and each ray in a generator system 
describing V satisfies all linear inequalities in a constraint system describing 
Q. As a further example, the time elapse operation specified in Section (H can 
be implemented using the generator systems for the argument polyhedra [70] . 
That is a generator system for the polyhedron V /* Q can be obtained by 
adopting the same set of generating points as V and by defining its set of 
rays as the union of the set of generating rays for V with the set of all the 
generators (both points and rays) for Q. 

As seen in Section [H in the context of the analysis of imperative languages one 
of the most frequent statements is variable assignment, where the expression 
assigned is safely approximated by an affine relation ip CI"xR". The (direct 
or inverse) image of an affine relation can be naively computed by embedding 
the input polyhedron V C R n into the space M 2n , intersecting it with the con- 
straints defining ip and finally projecting the result back on W 1 . However, due 
to the moves to/from a higher dimensional space, this approach suffers from 
significant overheads. Quite often, the expression assigned is a simple affine 
function of the variables' values and can thus be exactly modeled by comput- 
ing the image of a single-update affine function. With the double description 
method, the images of affine functions are much more efficiently computed 
by applying them directly to the generators of the argument polyhedron. A 
dual approach, using the constraint description of the polyhedron, allows for 
the computation of the preimages of affine functions, which can be of inter- 
est for a backward semantic construction, where the initial values of program 
variables are approximated starting from their final values. Similar efficiency 
arguments motivate the study of specific implementations for single-update 
bounded affine relations and other special subclasses of affine relations. 



7.2 Widening and Narrowing 

The first widening operator for the domain of convex polyhedra, the so-called 
standard widening proposed in [43] and refined in [65], can be informally de- 
scribed as follows: suppose that in the post-fixpoint iteration sequence we 
compute as successive iterates the polyhedra V% and Vi+i; then, the widen- 
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ing keeps all and only the constraints defining Vi that are also satisfied by 
Vi+i. This simple idea, which is basically borrowed from the widening oper- 
ator defined on the domain of intervals [37], is quite effective in ensuring the 
termination of the analysis (the number of constraints decreases at each iter- 
ation); by avoiding the application of the widening in the first few iterations 
of the analysis [35] and/or by applying the "widening up-to" technique of |66j, 
it also provides, in the main, an adequate level of precision. 

Some application fields, however, are particularly sensitive to the precision of 
the deduced numerical information, to the point that some authors propose to 
give up the termination guarantee and use so-called extrapolation operators: 
examples include the operators defined in [75] and [77], as well as the proposals 
in [29] and [46] for sets of polyhedra and the heuristics sketched in |26j. 

In |12| this precision problem is reconsidered in a more general context and a 
framework is proposed that is able to improve upon the precision of a given 
widening while keeping the termination guarantee. The approach, which builds 
on theoretical results put forward in work on termination analysis, combines an 
existing widening operator, whose termination guarantee should be formally 
certifiable, with an arbitrary number of precision improving heuristics. Its 
feasibility was demonstrated by instantiating the framework so as to produce a 
new widening on polyhedra improving upon the precision of [65] in a significant 
percentage of benchmarks. 

For the more challenging case of an abstract domain obtained by the finite 
powerset domain construction, several generic schemes of widenings have been 
proposed in [TTj that are able to "lift" a widening defined on the base-level do- 
main without compromising its termination guarantee. The instantiation of 
such a generic approach led to the definition of the first non-trivial and prov- 
ably correct widenings on a domain of finite sets of convex polyhedra. Being 
highly parametric, the widening schemes proposed in [T7] can be instantiated 
according to the needs of the specific application, as done in [63j. One of 
the heuristic approaches adopted in [17] to control the precision/complexity 
trade-off of the widenings, originally proposed in [29], attempts at reducing 
the cardinality of a polyhedral collection by merging two of its elements when- 
ever their set union happens to be a convex polyhedron. The implementation 
of such a heuristic could significantly benefit from the results and algorithms 
presented in [221124] . 

It is also worth mentioning that, once a post-fixpoint approximation has been 
obtained by means of an upward iteration sequence with widening, its precision 
can be improved by means of a downward iteration, possibly using a narrowing 
operator [37 38 40, 41]. To the best of our knowledge, no narrowing has ever 
been defined on the domain of convex polyhedra: applications simply stop the 
downward computation after a small number of iterations. 
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7.3 Not Necessarily Closed Convex Polyhedra 

Most static analysis applications computing linear inequality relations between 
program variables consider the domain CP n of topologically closed polyhedra. 
One of the underlying motivations is that sometimes (e.g., when working with 
integer valued variables only) strict inequalities can be filtered away by suitable 
syntactic manipulations; even when this is not the case, the topological closure 
approximation may be interpreted as a quick and practical workaround to the 
fact that some software libraries do not fully support computations on NNC 
polyhedra. However, there are applications [4 33 69 70] where the ability of 
encoding and propagating strict inequalities might be crucial for the usefulness 
of the final results. 

The first proposal for a systematic implementation of strict inequalities in 
a software library based on the double description method was put forward 
in [69]: a syntactic translation embeds an n-dimensional NNC polyhedron 
? GP„ into an (n + l)-dimensional closed polyhedron 7Z G CP n+ i, by adding 
a single slack variable e, satisfying the additional side constraints < e < 1. 
Namely, any strict inequality constraint (a, x) > b is translated into the non- 
strict inequality constraint (a, x) — e > b. The computation is thus performed 
on the closed representation 1Z G CP n+ i, with only minor adaptations to the 
basic algorithms so as to also take into account the implicit strict constraint 
e > 0. 

While this idea is quite effective, the resulting software library no longer en- 
joys all of the properties of the underlying double description implementation: 
NNC polyhedra cannot be suitably described using generator systems, and 
the geometric intuitions are lost under the "implementation details." These 
problems motivated the studies in [T3|14|19] . where a proper generalization of 
the double description method to NNC polyhedra was proposed. The main im- 
provement was the identification of the closure point as a new kind of generator 
for NNC polyhedra, leading to the following result generalizing Theorem 17.11 : 

Theorem 7.2 The set V C W 1 is an NNC polyhedron if and only if there 
exist finite sets R,P,C C M. n of cardinality r, p and c such that ^ R and 

V = \ Rp + Pit + C7 Gf 



The new condition tt ^ ensures that at least one of the points of P plays 
an active role in any convex combination of the vectors of P and C. As a 
consequence, the vectors of C are closure points of V, i.e., points that belong 
to the topological closure of V, but may not belong to V itself. 



p G M t ;,tt G R p + ,ir^ 0,7 G R° + , 
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Thanks to the introduction of (strict inequalities and) closure points, most of 
the pros of the double description method now also apply to the domain of 
NNC polyhedra: simpler, higher-level implementations of operations on NNC 
polyhedra can be specified, reasoned about and justified in terms of any one of 
the two dual descriptions; important implementation issues (such as the need 
to identify and remove all kinds of redundancies in the descriptions [H|T9]) 
can be provided with proper solutions; different lower-level encodings (e.g., an 
alternative management of the slack variable [T3~yT4] ) can be investigated and 
experimented with, without affecting the user of the software library. It would 
be interesting, from both a theoretical and practical point of view, to provide a 
more direct encoding of NNC polyhedra, i.e., one that is not based on the use 
of slack variables; this requires the specification and the corresponding proof 
of correctness of a direct NNC conversion algorithm, potentially achieving a 
major efficiency improvement. 



8 Conclusion 



In the field of automatic analysis and verification of software and hardware 
systems, approximate reasoning on numerical quantities is crucial. As first 
recognized in 1978 [43], polyhedral computation algorithms can be used for the 
automatic inference of numerical assertions that correctly (though usually not 
completely) characterize the behavior of a system at some level of abstraction. 

Until the end of the 1990's these techniques were not in widespread use, mainly 
due to the unavailability of robust and efficient implementations of convex 
polyhedra. As far as we know, the first published libraries of polyhedral al- 
gorithms suitable for analysis and verification purposes have been PolyliblH 
released in 1995, written by Wilde at IRIS A |113j and based on earlier work by 
Le Verge [85], and the polyhedra library of POLINE (POLyhedra INtegrated 
Environment) written by Halbwachs and Proy at Verimag and also released 
in 1995. Both libraries used machine integers to represent the coefficients of 
linear equalities and inequalities, something that could easily result into (unde- 
tected) overflows. While Polylib provided only a fraction of the functionalities 
offered by POLINE's library (which offered, among other things, support for 
NNC polyhedra), it was available in source format. The POLINE's library, in- 
stead, was distributed only in binary form for the Sun-4 platform (freely, until 
about the year 1996; under rather restrictive conditions afterward). POLINE 
included also a system called POLKA (POLyhedra desK cAlculator) and an 
analyzer for linear hybrid automata. A variation of a subset of POLINE's 
library was incorporated into the HyTech tool [76jjZI 



http : //www. ee .byu . edu/f aculty/wilde/polyhedra.html 



http : //embedded . eecs . berkeley . edu/research/hytech/ 
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The work of Wilde and Le Verge, which was extended by Loechner [87], led to 
the creation of PolyLib\^\ The New Polka library by Jeannet[l[] first released 
in 2000 and originally based on both IRISA's Polylib and POLINE's library, 
incorporates the idea — suggested by Fukuda and Prodon [57] — of lexico- 
graphically sorting the matrices representing constraints and generators. New 
Polka, which supports both closed and NNC polyhedra, together with Mine's 
Octagon Abstract Domain Library [90ll9lp1 and an interval library called 
ITV, is now included in the APRON library^ Finally, the Parma Polyhe- 
dra Library (PPL), initially inspired by New Polka and first released in 2001, 



is developed and maintained by the authors of this paperJiU The PPL sup- 
ports both closed and NNC polyhedra, bounding boxes, bounded difference 
and octagonal shapes, grids and combinations of the above including the finite 
powerset construction [T5"fT8] . 



The above libraries have all been designed specifically for applications of analy- 
sis and verification such as those described in this paper. However, two libraries 
that were designed for solving vertex enumeration/convex hull problems have 
successfully been used in static analysis and computer-aided verification tools: 



Fukuda's cddlib 
and Irslib 
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an implementation of the double description method [92J; 



the implementation by Avis of the reverse search algorithm [6]. 



All the libraries mentioned in the last two paragraphs are distributed under 
free software licenses and support the use of unbounded numeric coefficients. 
This, together with the ever increasing available computing power and the 
growing interest in ensuring the correctness of critical systems, has caused, in 
the 2000's, the continuous emergence of new tools and applications of poly- 
hedral computations in the area of formal methods. As a consequence, this is 
much more of a new beginning than an end to research in this area. As ex- 
plained in Sections [6] and [71 several open issues remain. Most of them have to 
do with the need for effectively managing the complexity-precision trade-off: 
the encouraging results obtained with today's tools are pushing us to apply 
them to more complex systems for a possibly more precise analysis and/or 
verification of more complex properties. 
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http : //icps . u- strasbg . f r/polylib/ 



http : //pop-art . inrialpes .f r/people/bjeannet/newpolka/index.html 
http : //www . di . ens . f r/~mine/oct/ 



http : //apron . cri . ensmp . f r/library/ 



http : //www. cs .unipr . it/ppl 



http : //www. if or .math. ethz . ch/~f ukuda/ cdd_home/ 



http : //cgm. cs .mcgill . ca/~avis/C/lrs .html 
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